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Abstract 

The spin-1 Ising (BEG) model with the nearest-neighbour bilinear and 
biquadratic interactions and single-ion anisotropy is simulated on a cellular 
automaton which improved from the Creutz cellular automaton(CCA) for 
a simple cubic lattice. The simulations have been made for several sets of 
parameters K/J and D/J in the -3 < D/J < and -1 < K/J < 
parameter regions. The re-entrant and double re-entrant phase transitions 
of the BEG model are determined from the temperature variations of the 
thermodynamic quantities (M, Q and x )• The phase diagrams characterizing 
phase transitions are compared with those obtained from other methods. 

Keywords: spin-1 Ising model; Creutz cellular automaton; re-entrant 
phase transition; simple cubic lattice. 
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1 . Introduction 

The Blume-Emery-Griffiths (BEG) modeP was originally introduced in 
order to explain the phase separation and superfluidity in the He'^-He^ mix- 
tures. Subsequently, the model was used in the description of a variety of dif- 
ferent physical phenomena such as multicomponent fluids^, microemulsions^, 
and semiconductors alloys'', etc. 

The Hamiltonian of the BEG model is given by, 

H,^jJ2s,S, + Kj2sfS] + Dj2s! (1) 

<ij> <ij> i 

where Sj = ±1, and < ij > denotes summation over all nearest-neighboring 
(nn) spin pairs on a simple cubic lattice. The parameters of J and K are 
the bilinear and biquadratic interaction energies, respectively and D is the 
single-ion anisotropy constant. The three-dimensional BEG model has been 
extensively studied by different techniques, using the mean-field approxima- 
tion (MFA)^^'^"'''^ effective- field theory (^"-'^^^ two-particle cluster approxima- 
tion (TPCA)^^, Bethe approximation^^, high-temperature scries expansion^"', 
renormalization group theory^^, Monte Carlo simulations*^^^'^^"^^-*, linear chain 
approximation^^^'^^) and cellular automaton^^^'^^^ 

In this paper we studied the three-dimensional BEG model using an im- 
proved heating algorithm from the Creutz Cellular Automaton (CCA) for 
simple cubic lattice. The CCA algorithm is a microcanonical algorithm 
interpolating between the canonical Monte Carlo and molecular dynamics 
techniques on a cellular automaton, and it was first introduced by Creutz^^. 
In the previous papers*^^°'^^'^^~^^\ the CCA algorithm and improved algo- 
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rithms from CCA were used to study the critical behavior of the different 
Ising model Hamiltonians on the two and three-dimensions. It was shown 
that they have successfully produced the critical behavior of the models. 

The BEG model has a complicated phase diagrams and has several kinds 
of phase transitions, such as re-entrant and double re-entrant transitions 
for some values of the model parameters on the three-dimensional lattices. 
However, there exists the differences between the phase diagrams in the 
{d — D/J, k — K/J, t — kBT/zJ) plane is obtained for certain parame- 
ter values^^'^^'^^'^^^ , where z is the coordination number. Recently, TPCA^^ 
calculations of the three dimensional BEG model with the ferromagnetic 
bilinear interaction show that the phase diagrams are qualitatively different 
from the phase diagrams obtained with the MFA^ and Bethe approximation^^ 
for the k — —0.5 parameter. In addition, the order parameters exhibit a dif- 
ferent critical behavior for the double re-entrant phase transition. Thus, the 
alternative results are needed to check the differences in results. The aim of 
the present paper is to study re-entrant behaviors of the BEG model in the 
case of J < and to obtain the phase diagrams for several sets of param- 
eters k and d in the —3 < o? < and — 1 < A; < parameter regions. For 
this purpose, the thermodynamic quantities are computed on a simple cubic 
lattice with linear dimension L = 18 and the finite critical temperatures are 
estimated from the maxima of the susceptibility. The model is explained in 
Sec. 2, the results are discussed in Sec. 3 and a conclusion is given in Sec. 4. 

2. Model 

Three variables are associated with each site of the lattice. The value of 
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each sites are determined from its value and those of its nearest-neighbors 
at the previous time step. The updating rule, which defines a deterministic 
cellular automaton, is as follows: Of the three variables on each site, the first 
one is Ising spin Bi. Its value may be or 1 or 2. The Ising spin energy for 
the model is given by Eq.l. In Eq.l, Si — Bi — 1. The second variable is for 
momentum variable conjugate to the spin (the demon). The kinetic energy 
associated with the demon, H^, is an integer, which equal to the changing in 
the Ising spin energy for the any spin flip and its values lie in the interval (0, 
m). The upper limit of the interval, m, is equal to 24 J. The total energy 

H^Hj + Hk (2) 

is conserved. 

The third variable provides a checkerboard style updating, and so it allows 
the simulation of the Ising model on a cellular automaton. The black sites 
of the checkerboard are updated and then their color is changed into white; 
white sites are changed into black without being updated. The updating 
rules for the spin and the momentum variables are as follows: For a site to 
be updated its spin is changed one of the other two states with 1/2 probability 
and the change is transferable to or from the momentum variable associated 
with this site, such that the total energy H is conserved, then this change is 
done and the momentum is appropriately changed. Otherwise the spin and 
the momentum are not changed. 

For a given total energy the system temperature is obtained from the 
average value of kinetic energy, which is given by: 
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m 

<E>=^ (3) 

Q-nJ/kT 

n=0 

where E — H^- The expectation value in Eq.3 is a average over the lattice 
and the number of the time steps. Because of the third variable, the algo- 
rithm requires two time steps to give every spin of the lattice a change to 
change. Thus, in comparison to ordinary Monte Carlo simulations, two steps 
correspond to one full sweep over the system variables. 

The heating algorithm is divided into two basic parts, initialization proce- 
dure and the taking of measurements. In the initialization procedure, firstly, 
all spins in the lattice sites take the ferromagnetic ordered structure (tt) and 
the kinetic energy per site which is equal to the maximum change in the Ising 
spin energy for the any spin flip is given to the certain percent of the lattice 
via the second variables. This configuration is run during the 10.000 cellular 
automaton time steps. At the end of the this step, the configuration in the 
ordered structure at the low temperature is obtained. In the next steps last 
configuration in the ordered structure has been chosen as a starting configu- 
ration for the heating run. Rather than resetting the starting configuration 
at each energy, it is used the final configuration at a given energy as the 
starting point for the next. During the heating cycle, energy is added to the 
system through the second variables (Hk) after 1.000.000 cellular automaton 
steps. 

3. Results and discussion 

The three-dimensional BEG model is simulated with the heating algo- 



5 



rithm which improved from the Creutz Cellular Automaton. The simulations 
were executed on simple cubic lattice LxLxL of the linear dimensions L = 18 
with periodic boundary conditions. The computed values of the quantities 
are averages over the lattice and over the number of time steps (1.000.000) 
with discard of the first 100.000 time steps during which the cellular au- 
tomaton develops. Hereafter we shall use the terminology of Ref. 13 for the 
definition of phases: F, the ferromagnetic phase (M 7^ 0, Q 7^ 2/3); Q, the 
quadrupolar phase (M = 0, (5^^2/3). 

In the ground state, the phase diagram of the BEG model on the (d, k) 
plane for — 1 < < parameter region is shown in Fig.l. The region of the 
perfect zero ordering and ferromagnetic ordering are separated by the line^^ 
k — —1 — d/S represented with the solid line. The types of phase transition 
(PT) obtained from the CCA calculations are also shown in Fig.l. In the 
ferromagnetic ordering region, the Q F phase transitions occur according 
to the values of k and d. However, the re-entrant Q ^ F ^ Q and double 
re-entrant Q ^ F ^ Q ^ F phase transitions take place for some values of 
k and d near the phase boundary. 

For the reproduce of the phase diagram using the CCA simulations, the 
temperature variation of the order parameters (M and Q) and the suscep- 
tibility (x) are calculated at d = 0, —0.25, —1 and —2.5 for several k values 
on the finite lattice with L — IS. For d — 0, the order parameters and the 
susceptibility data are shown in Fig. 2 (a) and (b). There is the second order 
phase transition from quadrupolar phase to ferromagnetic phase {Q — > F) at 
d = for A; > — 1. In this region, the susceptibility data shows characteristic 
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peaks at the critical temperature. However ior k — —1 the value of the mag- 
netization at the lowest temperature is much lower than the magnetization 
for the k > —1 values. 

The order parameters M and Q obtained at o? = —0.25 for different 
k values are shown in Fig.3 (a) and (b). For k > —0.90, there is only 
one second order Q ^ F phase transition in the system. In the interval 
—0.95 < k < —0.90, the re-entrant Q ^ F ^ Q phase transitions take 
place and both of them are second order. At the same time, the data of the 
susceptibility shows the two peaks which belong the Q ^ F and F ^ Q 
phase transitions in the —0.95 < k < —0.90 parameter region [Fig. 3(c)]. 

For d = —1, the temperature dependence of the order parameters and 
the susceptibility are given in Fig.4. As seen from Fig.4, while the second- 
order Q — > F phase transition occurs for k > —0.66, the double re-entrant 
Q ^ F ^ Q ^ F phase transitions occur for k = —0.66, —0.67 and 
—0.68. For k — —0.66, only the higher transition Q ^ F is second order 
and the others are first order. At the same time, the susceptibility for the 
lower and middle transitions show a sharp peaks at two different transition 
temperature [Fig.4 (c)]. Nevertheless, for k = —0.67 and —0.68 only the 
lower phase transition is first order because the order parameters exhibit a 
discontinuous behavior and the susceptibility exhibit a jump. According to 
these results, the middle phase transition F ^ Q in the double re-entrant 
Q ^ F ^ Q ^ F phase transition can be the first order and second order. 
This is in agreement with the results obtained by the MFA^ and TPCA^^ but 
disagreement with Bethe approximation(BA)^^ which gives only the lower 
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phase transition can be first order. 

In Fig. 5, the order parameters and the susceptibihty are shown for dif- 
ferent values of the k parameter at = —2.5. The re-entrant or double 
re-entrant PT's are not seen for d — —2.5. While the Q ^ F phase tran- 
sition is second order for k > —0.12, the Q ^ F phase transition is first 
order for k < —0.12. Nevertheless, the susceptibihty shows a characteris- 
tic behavior for k > —0.12 while a jump in the susceptibility appears for 
k < -0.12. 

Consequently, while the Q ^ F PT take place at d — ior k > —1, the 
re-entrant Q ^ F ^ Q and double re-entrant Q ^ F ^ Q ^ F PT's occur 
at c? = —0.25, —0.5, —1 and —1.5 for certain values of k. The representative 
phase diagrams relevant to the these phenomena are shown in Fig. 6. The 
phase diagrams obtained at d — and —1.5 represent similar behavior with 
the Bethe approximation results^^ illustrated in Fig. 6. 

Finally, we obtain the {d, t) phase diagram for k = —0.5 [Fig. 7]. It is quite 
different from the phase diagram obtained with different methods ^^'^^^ at 
k — —0.5. In the phase diagram obtained with the MFA^ and renormalization 
group theory^^ for k = —0.5, the second order PT line terminates at the 
critical end point E at the first order PT line and the first order PT line 
terminates at the critical point C inside the ferromagnetic phase. On the 
other hand, there is a tricritical point in the phase diagram obtained with 
TPCA^^ for k = —0.5. As seen Fig. 7, there is no critical point C in the 
ferromagnetic phase and only a tricritical point exists in the phase diagram 
with CCA calculations as in TPCA. 
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4. Conclusion 

The BEG model is simulated using the heating algorithm of the cellular 
automaton for simple cubic lattice. The CCA calculations show that the BEG 
model exhibits the re-entrant and double re-entrant phase transitions for 
some values of the k and d parameters near the phase boundary as expected. 
According to the Bethe approximation^^, the lower PT in the double re- 
entrant Q ^ F Q ^ F PT is of the first order and the other PT's are 
of the second order. Contrast to this, the MFA^ and TPCA^^ calculations 
suggest that the middle PT can be also the first order. The CCA results 
are compatible with the MFA and TPCA results at this point. We have also 
reconstructed the phase diagrams in the {d, k, t) parameters space using the 
CCA simulations. The phase diagrams in the (A;, t) plane for d = and —1.5 
are similar to the phase diagrams obtained by the BA^^ and TPCA^^ at those 
parameters. Although the MFA^ and renormalization group theory^^ predict 
that the {d, t) phase diagram for k — —0.5 has an end point at the first 
order phase transition line and a critical point in the ferromagnetic phase, 
the CCA calculations for k — —0.5 indicate that there is only a tricritical 
point as the results of TPCA^^. 
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Figure Captions 

Fig.l. The phase diagram (d, k plane) in the ground state and the types 
of phase transitions obtained on a simple cubic lattice. 

Fig. 2. The temperature dependence of (a) the order parameters M and 
Q, (b) the susceptibility x for d = 0. 

Fig. 3. The temperature dependence of (a) the order parameter M, b) Q, 
(c) the susceptibility x for (i = —0.25. 

Fig.4. The temperature dependence of (a) the order parameters M, b) 
Q, (c) the susceptibility x for d = — 1. 

Fig. 5. The temperature dependence of (a) the order parameters M and 
Q, (b) the susceptibility x for o? = —2.5. 

Fig. 6. The phase diagrams in the {k, t) plane for values of the c? = 0, 
-0.25, -0.5, -1 and -1.5. The closed and open symbols represent the second- 
order and first-order phase transitions for the CCA calculations, respectively. 
The dashed and thick dashed lines indicate the second-order and first-order 
phase transitions for the Bethe approximation (BA)^^ at the d — Q and —1.5 
parameter values. 
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Fig. 7. The phase diagram in the {d, t) plane for k — —0.5. The sohd hne 
and dotted hne represent the second-order and first-order phase transitions 
for TPCA^^ calculations. The closed and open squares represent the second- 
order and first-order phase transitions for CCA calculations, respectively. 
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